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We solve the coupled Wong Yang-Mills equations for both 17(1) and SU{2) gange groups and 
anisotropic particle momentum distributions numerically on a lattice. For weak fields with initial 
energy density much smaller than that of the particles we confirm the existence of plasma instabilities 
and of exponential growth of the fields which has been discussed previously. Also, the SU{2) case 
is qualitatively similar to 17(1), and we do find significant “abelianization” of the non-Abelian fields 
during the period of exponential growth. However, the effect nearly disappears when the fields are 
strong. This is because of the very rapid isotropization of the particle momenta by deflection in a 
strong field on time scales comparable to that for the development of Yang-Mills instabilities. This 
mechanism for isotropization may lead to smaller entropy increase than collisions and multiplication 
of hard gluons, which is interesting for the phenomenology of high-energy heavy-ion collisions. 

PACS numbers: 12.38.Mh,24.85.-|-p,25.75.-q,52.35.Qz 


High-energy heavy-ion collisions release a large amount 
of partons from the wavefunctions of the colliding nuclei. 
Partons with large transverse momenta originate from 
high-Q^ hard interactions which can be computed from 
perturbative QCD P|. On the other hand, partons with 
“small” transverse momenta on the order of the so-called 
saturation momentum Qs (given by the square root of 
the total color charge density per unit rapidity and unit 
area in the incoming nuclei) are much more abundant if 
Qs ^ Aqcd and are better viewed as a classical non- 
Abelian field 0. 

If the presence of the soft classical field is neglected, 
which amounts to assuming that Qs ~ Aqcd, the time- 
evolution of the hard partons after they come on-shell can 
be studied by means of the Boltzmann equation with a 
collision kernel, which is the so-called parton-cascade ap¬ 
proach iS The collision kernel could be truncated at 
the level of elastic binary collisions (perhaps with a sum¬ 
mation of time-like and space-like parton showers in the 
leading logarithmic approximation 0 ); recently, an at¬ 
tempt to fully include 2^3 processes beyond the relax¬ 
ation time and leading-logarithmic approximations has 
also been made [^. 

On the other hand, for large nuclei and at high ener¬ 
gies the saturation scale Qs is expected to grow much 
larger than Aqcd 00 and so the presence of the clas¬ 
sical field can no longer be neglected. The “bottom-up 
scenario” 0 generalizes the parton cascade description of 
the time-evolution after the collision to include the soft 
classical modes, too. Soft gluon radiation is found to be 
the dominant process leading to equilibration 000 
(see also papers by Wong in Q). 

Recently, it has been argued that collective processes 
due to the soft gauge field should be taken into account. 
Specifically, QCD plasma instabilities may develop due to 
anisotropic distributions of released hard parto ns ISl and 
modify the “bottom-up scenario” significantly [lj|. The 
hard loop effective action for anisotropic hard modes was 


formulated in [m and unstable soft modes were analyzed 
in Numerical studies of its static limit revealed 
the interesting tendency of the non-Abelian gauge fields 
to “abelianize” during the stage of instability in the sense 
that locally commutators become much smaller than the 
fields themselves (see below). The “abelianization” has 
also been seen in solutions of the full non-linear hard 
loop effective action 0 . It is argued that, because of 
abelianization, non-Abelian effects should not cause in¬ 
stabilities to saturate; rather, similarly to the Abelian 
case, the fields should continue to grow until their en¬ 
ergy density becomes comparable to that of the hard 
modes [lMllJ,ll3, i-e. until the growing fields begin to 
have a significant effect on the dynamics of the particles. 

It is interesting to note the following difference between 
isotropization by propagation of particles in a strong ran¬ 
dom field versus that via scattering and gluon multi¬ 
plication. Namely, in the absence of a collision kernel 
the entropy of any specific initial condition is conserved, 
while the standard parton cascade approach produces ad¬ 
ditional entropy 00 An ensemble average over suffi¬ 
ciently random initial field configurations can neverthe¬ 
less increase the entropy of the soft modes by a moderate 
(logarithmic) amount; this follows from the equivalence 
of the averaged classical field description to a Boltzmann 
equation to leading and subleading orders in the occupa¬ 
tion number 0 

In heavy-ion collisions, it might not be necessary to 
achieve “true” thermalization in the sense of maximizing 
the entropy during the first few fm/c of the reaction; 
isotropization could be sufficient In fact, data 0| 
from RHIC indicate that the number of charged particles 
per participant in the final state is only ~ 30% lower in 
central d-|-Au collisions than it is in central Au-|-Au. This 
perhaps indicates that the equilibration process expected 
to occur in Au-|-Au (but not in d-|-Au) does not produce 
a large amount of entropy |l 9|| . Hence, the mechanism of 
isotropization of particles via strong fields could be very 
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interesting for the phenomenology of heavy-ion collisions. 

In this letter we solve the classical transport equation 
for hard gluons with non-Abelian color charge (5“ in the 
collisionless approximation [iof . 

- gQ'^F^.d; - gfabcAlQ^dQa]f{x,p, Q) = 0 . (1) 
It is coupled to the Yang-Mills equation 

D^Ff^-=f = g J ^dQQv''f{x,p,Q), ( 2 ) 

where f{x,p, Q) denotes the one-particle phase space dis¬ 
tribution function |^. These equations were shown to 
reproduce the “hard thermal loop” effective action M 
near equilibrium. If fluctuations on top of the mean 
fields are not neglected, one obtains a collision term 
from their moments |2l| . The same set of transport 
equations were also derived within the worldline for¬ 
malism for the one loop effective action of QCD [ 2 ^ : 
the emergence of classical transport from a quantum ki¬ 
netic equation derived within the closed-time-path for¬ 
malism was discussed in ref. |2,*l|. For recent reviews 
of semi-classical transport theory for non-Abelian plas¬ 
mas see ref. |2^. Furthermore, we refer to ref. |2^ for a 
study of particle production and propagation in Abelian 
fields, including back-reaction and collisions in the re¬ 
laxation time approximation. The specific point of the 
present letter, however, is to study possible non-Abelian 
plasma instabilities due to anisotropic particle distribu¬ 
tions HESIilElElElEi 

We employ the test-particle method |2^ , replacing the 
continuous distribution f{x,p,Q) by a large number of 
test particles: 

f{x,p,Q) = -^'^S{x-Xi){2TrfS{p-p^)6(Qi-Q), 

i 

( 3 ) 

where and Pi are the coordinates of an individual test 
particle. This Ansatz leads to Wong’s equations [2(11^ 

^ = gQ-{E^ + v,xB^), (5) 

^ = igv'^[Af„Qi\ (6) 

for the i-th test particle [28l| . 

The time evolution of the Yang-Mills field can be fol¬ 
lowed by the standard Hamiltonian method . Numer¬ 
ical techniques to solve the classical field equations cou¬ 
pled to particles have been developed in ref. [s^]. Our 
update algorithm is closely related to the one explained 
there which generalizes the Abelian version of the charge 
conservation method in particle simulations [IJl- 

In the following, we assume that the fields only depend 
on time and on one spatial coordinate, x, which reduces 


the Yang-Mills equations to 1-1-1 dimensions. The hard 
modes represented by classical particles are allowed to 
propagate in three spatial dimensions. For simplicity, we 
also restrict ourselves to the case without expansion here; 
the more realistic case with longitudinal expansion 
will be addressed in the future. 

The initial anisotropic phase-space distribution of hard 
gluons is taken to be 

/(p, x) oc e-VPl+pl/p^^r^sip,) . (7) 

This represents a quasi-thermal distribution in two di¬ 
mensions, with “temperature” phard which now takes over 
the role of the saturation momentum mentioned above. 
We have checked explicitly that no instability occurs 
when the particle distribution is taken to be isotropic. 

The initial field amplitudes are sampled from a Gaus¬ 
sian distribution with a width tuned to a given initial 
energy density. We solve the Yang-Mills equations in 
A° = 0 gauge and also set A = 0 (i.e. all gauge links 
=1) at time t = 0; the initial electric field is taken to be 
polarized in a random direction transverse to the x-axis. 
Gauss’ law is then used to obtain the initial charge dis¬ 
tribution. All results shown below were obtained using a 
lattice with N = 512 sites; we have checked the numeri¬ 
cal accuracy by comparing to iV = 256, 1024 lattices (for 
the same physical parameters) and by monitoring conser¬ 
vation of the total energy and of Gauss’ law. The total 
energy was conserved to within 5 x 10“^ (5 x 10“^) over 
the course of the simulations for the weak (strong) field 
initial conditions, and the maximal violation of Gauss’ 
law was 10“® for SU{2) and 10“^^ for U{1) (in lattice 
units). 

Before coming to our results, we also comment on the 
occurence of “anomalous Gherenkov radiation”. This 
corresponds to anomalous hard radiation from soft modes 
which may occur for simulations on a discrete lattice, 
as the dispersion relation of the fields may contain real 
space-like modes. For example, taking uj(k) = k for 
free fields in the continuum leads to the dispersion re¬ 
lation u;{k) = 2| sin(fca/2)|/a on a one-dimensional lat¬ 
tice (|fc| < n/a). Gonsequently, hard field modes with 
fc ~ 1/a would then get populated on the lattice be¬ 
cause their “mass” is imaginary. The situation 

could perhaps be improved by employing higher-order 
discretization schemes for the Yang-Mills action or by 
damping hard modes exponentially at t = 0. However, 
we have not done so at present. Our tests with different 
lattices do not indicate a significant dependence of either 
the growth or the saturation of the instability on the 
lattice spacing. Also, our solutions for isotropic particle 
momentum distributions appear to be stable when the 
number of test particles is taken to infinity. A possible 
physical reason for this observation could be that interac¬ 
tions among soft modes with fc <C 1/a and with large oc¬ 
cupation numbers, and interactions of those modes with 
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the particles, which give the largest contribution to the 
total energy (see below), actually dominate. Neverthe¬ 
less, this effect may deserve a more careful numerical 
study in the future. 

We first show results for a large separation of initial 
particle and field energy densities which should qualita¬ 
tively resemble the conditions studied in The 

results shown in Fig. ^ corresponds to a lattice of phys¬ 
ical size i = 40 fm, a hard scale Phard = 10 GeV and a 
particle density oinjg^ = 10 fm“^. 



FIG. 1: Time evolution of the kinetic (particle), magnetic 
and electric energy densities in GeV/fm® for U{1) and SU{2) 
gauge group, respectively. 

For the [/(I) case we observe a rapid exponential 
growth of the magnetic field energy density starting 
at about t/L « 0.1, turning into a slower growth at 
tjL « 0.5; at this point the magnetic fields have grown 
sufficiently to affect the particles which visibly start loos¬ 
ing energy. The electric field grows less rapidly and 
equipartitioning is not achieved within the depicted time 
interval. This indicates that the field strengths are still 
too high for linear response to apply. In the non-Abelian 
case the growth of the magnetic field saturates earlier, 
and the electric field has comparable strength by the end 
of the simulation. Also, it appears that the saturation of 
the magnetic instability occurs before it has a noticeable 
effect on the particles since their energy density is nearly 
constant. Nevertheless, at a purely qualitative level the 
C/(l) and SU(2] simulations are not vastly different, as 
anticipated in [iSj. 

This is analyzed further in Fig. [3 showing the growth 



FIG. 2: The average amplitude 4>rnis (in units of GeV) and the 
relative size C of commutators as a function of time; physical 
parameters as in Fig. 0 


of the rms average 


dx 

J _(A“A“ + A:a:) 


1/2 


( 8 ) 


and the average of the relative size of the field commuta¬ 
tor defined by 0 

dx \/Tr {{i[Ay,A^])'^) 

Jo L Tt{AI+AI) ■ 

The behavior of (prms is similar to that of the field energy 
density shown above. Initially, C is constant but then 
starts dropping exponentially when the magnetic insta¬ 
bility sets i^ indicating the partial “abelianization” of 
the fields [T^ H^ . The rate by which C drops in the 
intermediate stage is roughly comparable to the growth 
rate of (farms', also, the abelianization appears to stop after 
C dropped by about one order of magnitude, at about 
the same time when the exponential growth of the fields 
and of (/frms saturates. 

Finally, in Fig. |3 we show the time evolution of the 
longitudinal and transverse components of the energy- 
momentum tensor of the particles, i.e. the kinetic pres¬ 
sure. For both t/(l) and SU{2) we observe a rapid 
growth of the longitudinal pressure, which is zero ini¬ 
tially. Again, the rate is somewhat smaller for the non- 
Abelian case. The approach to “isotropization” of the 
kinetic pressure is clearly correlated to the stage of expo¬ 
nential growth of the soft fields seen before in . However, 
for both cases T^x remains significantly smaller than the 
transverse component for times ^ L. 

The initial conditions above were chosen such as to 
verify qualitatively the picture emerging in the hard loop 
approximation, where the field energy density is (and re¬ 
mains) much smaller than that of the par ticles and so 
the back-reaction can be neglected . In the 
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FIG. 3: Transverse and longitudinal components of the 
energy-momentum tensor of the particles for the simulation 
corresponding to Figs. Ill2l (weak field). 

color glass condensate model of high-energy collisions one 
does not expect such a strong separation of energy densi¬ 
ties, however. Since our numerical solution includes the 
back-reaction of the fields on the particles, we study the 
situation with stronger helds next. 

Specifically, the simulations below were performed 
with the following set of physical parameters: length 
L = 10 fm, hard scale phard = 1 GeV, particle density 
njg^ = 500 fm“^ and an initial field energy density of 
about 20 GeV/fm^. 



FIG. 4: Time evolution of the kinetic and field energy den¬ 
sities for strong initial fields. 


Fig. a shows the time evolution of the energy densi¬ 
ties for these initial conditions. This case clearly differs 
from the weak-field limit shown before. Over the time 
interval shown, the electric field energy density is prac¬ 
tically constant for both 17(1) and SU{2). The Abelian 
magnetic field does exhibit a slow growth, draining some 
energy from the particle reservoir. For SU{2), however, 
after a short initial growth the magnetic field energy de¬ 
creases again, to saturate pretty much at its initial value. 
Therefore, the particle energy density is also more or less 
constant over the depicted time interval. 



FIG. 5: Time evolution of (fiims and C in the strong field case. 

Fig. m confirms this observation via the ^rms observ¬ 
able: the initial growth saturates much earlier than be¬ 
fore. Similarly, the average commutator C stays constant 
for some time then drops by about a factor 5 (during the 
period where the magentic field drops !) and saturates 
at « 10%, which is an order of magnitude larger than for 
the weak field case from Fig. |3 This indicates a much 
smaller degree of “abelianization” for strong helds. 

Perhaps surprisingly. Fig. El nevertheless shows a very 
rapid isotropization of the kinetic pressure for both [/(I) 
and SU{2) (note that t/L = 0.1 corresponds to t=l fm 
in physical units for this lattice). Moreover, the degree 
of isotropization is much higher, i.e. the transverse and 
longitudinal pressures are closer than in Fig. El The very 
fast and nearly complete isotropization is, of course, the 
reason why held instabilities can not be sustained over 
a signihcant period of time in this case. It is caused by 
the bending of the particle trajectories in the strong held 
which is very different from conventional parton cascade 
transport with small-angle perturbative scattering (and 
no held). The random initial helds then cause a rapid 
isotropization of the particle momenta via eq. Note 
that this does not require hard modes in the helds, which 
indeed would violate the assumed separation of momen¬ 
tum scales, but large held amplitudes. 

In summary, we have studied instabilities in the cou¬ 
pled Wong Yang-Mills equations for strongly anisotropic 
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FIG. 6; Transverse and longitudinal components of the 
energy-momentum tensor of the particles for the simulation 
corresponding to Figs. 14151 (strong field). 

initial particle momentum distributions. For both 17(1) 
and SU{2) gauge groups we do observe a period of ex¬ 
ponential growth of the fields when their initial energy 
density is far less than that of the hard modes (parti¬ 
cles). This, in turn, leads to partial isotropization of the 
particle momentum distributions and of the kinetic pres¬ 
sure. Although we find somewhat smaller field growth 
and isotropization rates for the non-Abelian case, we 
nevertheless qualitatively confirm the picture developed 
in [HEi in that the non-Abelian fields “abelianize” ef¬ 
ficiently during the period of exponential growth. 

For large initial field amplitudes, corresponding to a 
smaller ratio of initial particle to field energy densities, 
our results are qualitatively different. We observe a very 
rapid isotropization of the particle momentum distribu¬ 
tions which is due to bending of their trajectories in 
the strong fields on a time-scale that is relevant for the 
physics of high-energy collisions. This, however, inhibits 
the development of instabilities of the Yang-Mills fields. 
Nevertheless, these results, too, suggest that the presence 
of the strong non-Abelian fields should be taken into ac¬ 
count to understand the process of isotropization in the 
early stages of high-energy collisions. 

We thank F. Gelis, M. Gyulassy, A. Mueller, R. Pis- 
arski and M. Strickland for helpful discussions and 
J. Lenaghan and M. Strickland for drawing our atten¬ 
tion to this subject. 

Note added: After this manuscript was submitted 
for publication, a paper appeared |3^ which presents an 
analytical discussion of a possible effective potential for 
anisotropic QGD plasmas beyond the hard loop approxi¬ 
mation. Also, a modified “bottom-up” scenario for gluon 
thermalization in high-energy heavy-ion collisions M 
and 3 -I- Id simulations of the full hard-loop effective the¬ 
ory appeared. 


[1] J. P. Blaizot and A. H. Mueller, Nucl. Phys. B 289, 847 
(1987); K. J. Eskola, K. Kajantie and J. Lindfors, Nucl. 
Phys. B 323, 37 (1989); X. N. Wang and M. Gyulassy, 
Phys. Rev. D 44, 3501 (1991). 

[2] L. McLerran and R. Venugopalan, Phys. Rev. D 49, 2233 
(1994); ibid. 49, 3352 (1994); Y. V. Kovchegov, ibid. 54, 
5463 (1996); ibid. 55, 5445 (1997); 

[3] K. Geiger and B. Muller, Nucl. Phys. B 369, 600 (1992); 
S. A. Bass, B. Muller and D. K. Srivastava, Phys. Lett. 
B 551, 277 (2003). 

[4] M. Gyulassy, Y. Pang and B. Zhang, Nucl. Phys. A 626, 
999 (1997); D. Molnar and M. Gyulassy, Phys. Rev. C 62, 
054907 (2000); J. Bjoraker and R. Venugopalan, Phys. 
Rev. C 63, 024609 (2001); S. Cheng et al, Phys. Rev. C 
65, 024901 (2002); Y. Nara, S. E. Vance and P. Csizma- 
dia, Phys. Lett. B 531, 209 (2002); 

for approaches based on the relaxation time approxima¬ 
tion to the Boltzmann equation see T. S. Biro, E. van 
Doom, B. Muller, M. H. Thoma and X. N. Wang, Phys. 
Rev. C 48, 1275 (1993); S. M. H. Wong, Nucl. Phys. A 
607, 442 (1996); Phys. Rev. C 56, 1075 (1997); Nucl. 
Phys. A 638, 527C (1998); G. C. Nayak, A. Dumitru, 
L. D. McLerran and W. Greiner, Nucl. Phys. A 687, 457 
(2001); A. Dumitru and M. Gyulassy, Phys. Lett. B 494, 
215 (2000); J. Serreau and D. Schiff, JHEP 0111, 039 
( 2001 ). 

[5] Z. Xu and C. Greiner, arXiv:hep-ph/0406278 

[6] A. H. Mueller, Nucl. Phys. B 558, 285 (1999). 

[7] R. Baier, A. H. Mueller, D. Schiff, and D. T. Son, Phys. 
Lett. B 502, 51 (2001). 

[8] G. R. Shin and B. Muller, J. Phys. G 28, 2643 (2002). 

[9] S. Mrowczynski, Phys. Lett. B 214, 587 (1988); Phys. 
Lett. B 314, 118 (1993); Phys. Rev. C 49, 2191 (1994); 
J. Randrup and S. Mrowczynski, Phys. Rev. C 68, 034909 
(2003). 

[10] P. Arnold, J. Lenaghan, and G. D. Moore, JHEP 0308, 
002 (2003). 

[11] S. Mrowczynski, A. Rebhan and M. Strickland, Phys. 
Rev. D 70, 025004 (2004). 

[12] P. Romatschke and M. Strickland, Phys. Rev. D 68, 
036004 (2003); Phys. Rev. D 69, 065005 (2004); 

[13] P. Arnold and J. Lenaghan, Phys. Rev. D 70, 114007 
(2004). 

[14] A. Rebhan, P. Romatschke and M. Strickland, 
arXiv:hep-ph/0412016 

[15] P. Arnold, J. Lenaghan, G. D. Moore and L. G. Yaffe, 
Phys. Rev. Lett. 94, 072302 (2005). 

[16] K. Geiger, Phys. Rev. D 46, 4986 (1992); S. M. H. Wong, 
Phys. Rev. C 54, 2588 (1996). 

[17] A. H. Mueller and D. T. Son, Phys. Lett. B 582, 279 
(2004). 

[18] see, for example. Pig. 3 in B. B. Back et al. [PHOBOS 
Collaboration], arXiv:nucl-ex/0409021 

[19] See also discussion following eq. (1) of S. Pratt, Nucl. 
Phys. A 715, 389 (2003). 

[20] P. F. Kelly, Q. Liu, C. Lucchesi and C. Manuel, Phys. 
Rev. Lett. 72, 3461 (1994); Phys. Rev. D 50, 4209 (1994); 
J. P. Blaizot and E. lancu, Nucl. Phys. B 557, 183 (1999). 

[21] D. F. Litim and C. Manuel, Phys. Rev. Lett. 82, 4981 
(1999); Nucl. Phys. B 562, 237 (1999). 

[22] J. Jalilian-Marian, S. Jeon, R. Venugopalan and 












6 


J. Wirstam, Phys. Rev. D 62, 045020 (2000). 

[23] Q. Wang, K. Redlich, H. Stocker and W. Greiner, Phys. 
Rev. Lett. 88, 132303 (2002). 

[24] J. P. Blaizot and E. lancu, Phys. Rept. 359, 355 (2002); 
D. F. Litim and C. Manuel, Phys. Rept. 364, 451 (2002). 

[25] D. V. Vinnik, A. V. Prozorkevich, S. A. Smolyan- 
sky, V. D. Toneev, M. B. Hecht, C. D. Roberts and 
S. M. Schmidt, Eur. Phys. J. C 22, 341 (2001). 

[26] G. F. Bertsch and S. Das Gupta, Phys. Rept. 160 (1988) 
189. 

[27] S. K. Wong, Nuovo Cim. A 65, 689 (1970); H. T. Elze 
and U. W. Heinz, Phys. Rept. 183, 81 (1989). 

[28] In our numerical calculations we redefine the Hamilto¬ 
nian H' — g^H, the fields = gA^ and the particle 
phase space density /' = g^ f in order to remove the 
gauge coupling from the classical theory. 

[29] J. Ambjorn, T. Askgaard, H. Porter and M. E. Shaposh- 
nikov, Nucl. Phys. B 353, 346 (1991). 

[30] G. R. Hu and B. Muller, Phys. Lett. B 409, 377 (1997); 


G. D. Moore, C. R. Hu and B. Muller, Phys. Rev. D 58, 
045001 (1998). 

[31] J. Villasenor and O. Buneman, Comp. Phys. Comm. 69, 
306 (1992). 

[32] A. Krasnitz and R. Venugopalan, Nucl. Phys. B557, 237 
(1999); Phys. Rev. Lett. 84, 4309 (2000); Phys. Rev. 
Lett. 86, 1717 (2001); A. Krasnitz, Y. Nara, and R. Venu¬ 
gopalan, Phys. Rev. Lett. 87, 192302 (2001); Nucl. Phys. 
A717, 268 (2003); Nucl. Phys. A 727, 427 (2003); Phys. 
Lett. B 554, 21 (2003); T. Lappi, Phys. Rev. C 67, 
054903 (2003). 

[33] C. Manuel and S. Mrowczynski, hep-ph/0504156 

[34] A. H. Mueller, A. 1. Shoshi and S. M. H. Wong, 
hep-ph/0505164 

[35] P. Arnold, G. D. Moore and L. G. Yaffe, 
arXiv:hep-ph/0505212 A. Rebhan, P. Romatschke 
and M. Strickland, arXiv:hep-ph/0505261 




FIG. 7: ADDITIONAL FIGURE: Time evolution of the field energy density on two different lattices (with the same physical 
size L); Np denotes the number of test particles per lattice site. 



FIG. 8: ADDITIONAL FIGURE: Same for the commutator. 










FIG. 9: ADDITIONAL FIGURE: time evolution of the field energy density for isotropic particle momentum distributions. 




